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ABSTRACT 


A  quantitative  analysis  of  the  Cut -Bar  method  of  measuring  the  thermal  con¬ 
ductivity  of  solids  is  performed.  The  mathematical  model,  which  corrects  for  the 
difference  in  heat  flux  in  the  specimen  and  reference  standard,  is  that  of  the  two 
dimensional  steady  heat  conduction  equation  applied  to  an  annulus  of  insulation.  The 
solution  is  presented  in  detail  and  found  to  be  comprised  of  two  physically  distinct  parts, 
a  conductivity  factor  and  a  geometrical  factor.  A  number  of  charts  and  graphs  are 
presented  for  clarification  as  to  the  nature  and  magnitude  the  relative  sizes  of  the 
various  components  will  have  on  the  accuracy  over  different  conductivity  ranges.  The 
complexity  of  the  geometrical  factor  required  a  digital  computer  programs  which 
are  included. 

Reference  is  made  to  a  similar  study,  performed  by  researchers  at  the 
National  Bureau  of  Standards.  It  is  found  that  the  differences  in  the  guard  tempera¬ 
ture  distribution  results  in  a  substantial  change  in  the  geometrical  factor. 
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I.  INTRODUCTION 


The  "Cut-Bar"  technique  of  measuring  thermal  conductivity  is  a  steady  state 
comparative  method  which  is  most  accurate  in  the  low  and  intermediate  temperature 
ranges  for  highly  conductive  materials.  The  schematic  of  the  elements  of  a  typical 
apparatus  is  illustrated  in  Figure  1.  It  consists  of  a  pair  of  meter  bars  or  discs 
between  which  a  specimen  is  interposed.  This  composite  bar  is  surrounded  by  an 
annulus  of  powder  insulation  which,  in  turn,  is  encased  in  a  heated  metal  cylinder 
acting  as  a  guard  against  extraneous  heat  losses .  An  axial  heat  flux  is  established 
by  an  isothermal  source  and  sink  at  either  end  of  the  assembly. 

For  a  first  approximation,  the  determination  of  the  specimen’s  conductivity  is 
rather  straightforward.  Fourier’s  equation  is  written  once  for  the  average  of 
the  meter  bar  valves  and  once  for  the  specimen: 


q/a|  =  k  s 

^  '  m  mm 


Q/A  =  KsSs 


(1) 

(2) 


where  Q/A  is  the  heat  flux,  k  the  thermal  conductivity,  and  S  the  temperature 
gradient  in  the  axial  direction.  By  assuming  the  heat  flow  and  cross-sectional 
areas  constant  throughout  the  composite  bar  the  following  relationship  exists: 

Ks  =  Km  [Sm/Ssl  :  (3) 

where  the  temperature  gradients  are  measured  and  the  conductivity  of  the  meter 
bar  is  known  either  by  using  a  known  standard  material  for  the  meter  bar  or  by 
previous  calibration  with  known  standards  as  specimens . 

Actually,  the  constant  heat  flow  assumption  is  in  error  for  two, reasons;  the 
radial  heat  exchange  with  the  guard  and  the  axial  shunting  exchange  with  the  powder 
insulation  caused  by  the  difference  in  conductivities  between  the  meter  bars  and 
specimen  (Figure  2).  A  more  realistic  expression  of  the  relative  heat  fluxes  is: 

KmSm  (1-y  m)  =  KgSs  (1-y  s)  (4) 


1 


where  y  is  the  ratio  of  the  heat  flow  crossing  the  interface  of  the  bar  and  insulation 
to  the  longitudinal  heat  flow  in  the  bar  if  there  were  no  losses  (value  in  the  bar 
a  differential  distance  from  the  source).  The  evaluation  of  y  for  a  particular  set  of 
conditions  may  be  done  either  analytically  or  experimentally.  This  report  is  con¬ 
cerned  with  an  analytical  solution. 
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FIGURE  1. 
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CUT-BAR  APPARATUS  WITH  A  LINEAR  GUARD 
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II.  ANALYSIS 


A  theoritical  evaluation  of  y  may  be  obtained  by  determining  the  heat  flow  across 
the  interface  of  the  bar  and  powder  insulation.  The  heat  conduction  equation  and 
boundary  conditions  for  the  annulus  of  insulation  are  assumed  to  be: 


p  =  a 
p  =  b 

(1  =  0 

M  =  1 


a2e 

ap5 


+  3_e_  +  s  2q 


p  3  p 


3  42 


=  G(z)  -  2 
Dl  -  U2 

_  H(z)  -  v  2 
U  !  -  U2 


(5) 


(6) 


Where  the  notation  may  be  realized  by  reference  to  Figure  3  and  the  nomenclature  section. 
The  general  solution  to  equation  (5)  is  of  course  known  (Reference  [1])  and  when  applied  to 
the  particular  boundary  conditions  the  temperature  distribution  throughout  the  insulation 
is  obtained.  This  distribution  expression  may  be  differentiated  to  obtain  the  radial 
gradient  which  in  turn  may  be  evaluated  at  the  A  radius  boundary  and  then  integrated 
along  the  A  axis  to  get  the  total  heat  flow  crossing  the  bar-insulation  interface. 

Dividing  this  result  by  the  longitudinal  heat  flow  in  the  bar  (assuming  no  losses) 
results  in  the  following  expression: 
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where  the  complete  details  of  the  above  mentioned  steps  are  presented  in  Appendix  A. 
Equation  (7)  may  be  seperated  into  two  terms  of , independent  physical  meaning,  the 
conductivity  and  geometrical  factors: 


Because  of  its  complexity  equation  (9)  has  been  programmed  and  is  usually  found  to 
converge  on  a  sufficiently  accurate  value  in  less  than  150  terms.  The  flow  chart  and 
programs  are  presented  in  Appendix  B.  For  a  given  apparatus  F  is  seen  to  be  a 
function  of  z  only.  However  since  the  specimen  length  (L)  effects  the  overall  length 
(w)  the  actual  values  of  F  (z)  must  be  determined  for  each  test.  In  Figure  4  a  curve 
for  particular  inputs  of  w  =  9.5  inches,  A  =  1  inch,  B  =  3  inches,  L  =  2  inches, 
represents  the  variation  of  the  geometrical  factor  (Fg)  with  dimensionless  distance  (y  ). 

,  A  reproduction  of  the  computer  program  output  can  be  seen  in  Figure  5  for  one- 
half  of  the  cut-bar.  Since  the  ends  of  the  cut -bar  and  guard  are  matched,  the 
radial  losses  and  gains  are  symmetrical  about  the  center  of  the  apparatus;  the 
geometrical  factor  F  (Z)  plotted  over  one -half  of  the  cut -bar  is  a  mirror  image  for 
the  other  half. 

Finding  Ft  for  the  apparatus  from  equation  (8)  and  using  F  (Z)  from  the  com- 
K  g 

puter,  the  same  curve  with  different  ordinate  scale  will  represent  y  (Z). 
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In  an  earlier  work  (Reference  2),  Flynn  carried  out  an  analysis  for  the  case  of  the 
outer  boundary  temperature  distribution  equal  to  the  inner  boundary;  that  is,  H(Z)  =  G(Z) 
A  comparison  with  this  ideal  case  is  made  in  Figure  6.  The  cut -bar  apparatus  with  a 
matched  guard  has  smaller  geometrical  factors  compared  with  a  similar  apparatus  and 
a  linear  guard  matched  at  the  ends,  (Figure  6).  The  two  cut -bar  apparatus  have  the 


same  dimensionless  ratios, 


L  _ 


=  3 


and  -Li-  =  4,  and  differ  only  in  the  thermal 
w  A 

guarding.  The  relative  flatness  of  the  matched  guard’s  correction  curves,  has  the 

advantage  of  not  having  to  average  the  correction  factor  between  thermocouples 

B 

discussed  in  the  example  in  Chapter  III.  Furthermore,  as  the  rat io-^— decreases  in 
the  linear  guard  apparatus,  the  geometrical  factor  will  increase  to  a  very  large  number. 
Conversely,  letting -r-de crease  in  the  matched  guard  apparatus,  the  geometrical  factor 
will  decrease,  and  in  the  limit  as-^-  approaches  one,  F^  approaches  zero.  There  phenomena 
are  better  illustrated  in  Figure  7,  which  is  a  plot  of  the  maximum  geometrical  factors  (at 
the  mid-horizontal  plane)  against  the  ratio  of  the  guard  and  bar  radii.  The  linear  guard’s 
increasing  geometrical  factor  with  decrease  of  may  be  explained  with  the  decreasing 
insulating  resistance,  the  radial  losses  begin  to  approach  the  axial  flux.  The  matched 

1  decreases  the  shunting 
loss  by  decreasing  the  medium  through  which  it  takes  place.  The  conditions  of  the  zero 
limit  are,  of  course,  quite  impractical. 

It  may  be  concluded  that  although  the  errors  associated  with  the  linear  guard  are  of 
greater  magnitude  and  variation  along  the  longitudinal  axis,  than  that  for  a  matched  guard, 
they  are  not  excessive  for  large  JL  ratios.  Consequently,  this  method  may  still  be  the 


guard  apparatus  has  a  minimum  of  radial  losses;  thus,  letting 


A 

preferred  technique;  particularly,  when  considering  the  additional  complexity  and  cost 
of  the  matched  guard’s  experimental  apparatus .  The  details  of  how  to  design  and 
correct  a  linear  guard  apparatus  are  discussed  in  Chapter  III. 
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FIGURE  5 
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COMPUTER  OUTPUT 


RESULTS  OF  POWER  LOSS  OVER  LENGTH  OF  METER- BAR 
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FIGURE  6 

COMPARISON  OF  THE  GEOMETRICAL  FACTORS  FOR  A  CUT-BAR 
APPARATUS  HAVING  A  MATCHED  GUARD  WITH  AN  APPARATUS 
HAVING  A  LINEAR  GUARD  MATCHED  AT  THE  ENDS 

The  dimensionless  variables  of  the 
apparatus  are:  L  =  3,  L  =  4 

W  A 


L„  G.  -  Linear  guard  matched  at  the  ends 
M.  G.  -  Matched  Guard 
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FIGURE  7 


GEOMETRICAL  FACTOR  AS  A  FUNCTION  OF  THE  RATIO 
OF  APPARATUS  DIAMETERS 


Plot  of  the  geometrical  factor  as  a  function  of 
the  ratio  of  guard  and  specimen  diameters  for  a 
matched  guard  and  a  linear  guard  matched  at  ends 
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III.  DESIGN 


The  computer  results  of  equation  (9)  are  also  plotted  as  design  charts .  The 
design  charts  can  be  used  to  find  optimum  dimensions  for  a  test  setup.  A  larger  than 
practical  range  of  variables  are  shown,  so  that  the  designer  can  be  sure  of  an  optimum. 
It  should  be  noted  that  in  the  development  of  the  design  charts,  the  only  dimensionless 
parameters  are  those  with  physical  implications . 

The  design  chart  (Figure  8)  plots  the  geometrical  factor  Fg  against  the  radii 

T)  L 

ratio  of  guard  to  bar-^-  ,  using  specimen  length  to  specimen  radius  and  overall 

length  i-,  as  parameters.  The  geometrical  factor  F  ,  should  be  kept  to  a  minimum 

®  B 

to  obtain  the  optimum  design.  It  is  evident  from  equation  (9)  that  as  -r-  approaches 

L  L  3 

1,  F  approaches  infinity  as  a  limit  for  any  -r-  and  .  However,  as  -j-  increases 

F^  rapidly  settles  down  to  various  asymptotic  values.  An  arbitrary  curve  designated 

®  L  L 

Bc  is  drawn  from  a  locus  of  points  in  the  family  of  curves  for  and  where 

the  slope  of  these  curves  have  an  absolute  value  of  0.01.  An  increase  in  the  ratio 

Bc 


passed  the  ratio  Bc 

a 


does  not  significantly  decrease  Fg.  Thus, 


A 


may  be 


considered  as  points  of  diminishing  return  and  are  the  recommended  value  for 
optimum  design. 

Figure  9  is  a  plot  of  F^  against  for  an  family  of  curves.  It  also  has  a 

®  t*  ^  B 

scale  on  the  right  ordinate  of  -^-c  .  If  the  designer  uses  ~-c  as  the  sizing  of  the 
annulus,  he  may  immediately  find  the  minimum  geometrical  factor  Fg,  for  whatever 
values  of — and  he  desires  .  It  may  be  seerftn  Figure  9,  that  as  approaches 
the  limit  — —  =  1,  F  -*  0  for  all  values  of  When  =  1,  the  cut -bar  is  a  solid 
bar  composed  of  the  specimen  material;  similarly,  when  L  =  0,  the  cut -bar  is  a  solid 
bar,  but  composed  of  meter  bar  material.  In  both  cases,  there  is  an  axial  temperature 
gradient  in  the  cut-bar  which  matches  the  axial  temperature  gradient  in  the  guard. 
Therefore,  F  must  be  zero,  because  there  is  neither  a  radial  heat  flow  nor  a  shunting 
flow  in  the  annulus . 

The  effect  of  varying  the  specimen  thermal  conductivity  on  the  thermal  con¬ 
ductivity  factor,  F^,  (Figure  10)  can  be  seen  by  using  equation  (8)  to  plot  F^  against 


K, 

Ki 


K 


for  various  ratios  of 
K  s 


m 


Ki 


The  plotted  quantities  go  through  zero  at 


K, 


Ki 


K 


m 


Ki 


then  approach  “  as  becomes  small,  and  asymptotically  approach  - 


K  becomes  large. 


Km 


as 


K. 

i 
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Since  a  designer  would  choose  a  geometrical  factor  which  would  be  small,  he  could 
multiply  the  ordinary  scale  of  Figure  10  to  obtain  y  .  The  ordinate  would  be  the 
measure  of  thermal  conductivity  of  materials  in  terms  of  fractional  power  gained  or 
lost  in  the  cut -bar  apparatus. 

If  it  is  necessary  to  measure  a  large  range  of  specimen  thermal  conductivities 
with  a  single  pair  of  meter-bars,  the  conductivity  of  the  meter-bar  should  be  in  the 
lower  end  of  the  range.  To  have  a  minimal  thermal  conductivity  factor  range,  for  a 
given  range  of  specimen  materials  using  only  a  single  pair  of  meter-bars,  it  is 
necessary  to  find  the  arithmetic  average  of  the  resistance  of  the  specimen  range. 

The  arithmetic  average  of  the  maximum  and  minimum  of  the  specimen  resistance 
range  will  be  used  as  the  optimum  resistance  for  the  meter -bar. 

Therefore, 


R 


m 


1/2 


s ^  max  + 


min 


(11) 


Substituting  the  equivalent  thermal  conductivity  for  resistance,  that  is,  R 


into  equation  (11) 
it  becomes, 


K 


m 


(k>)  (K's) 

'  7  max  \  / 


min 


/K\  +  /  K  \ 

s  1  max  I  s  1 


min 


(12) 


(Ko)  and  (KG)  are  the  maximum  and  minimum  specimen  thermal  con- 
'  max  v  a  mm 

ductivities  to  be  measured. 

The  maximum  fractional  power  is  approximately, 

1  1 


K. 

y  max  =  +  — ^ — 


(K 


s'  min 


max 


g 


(13) 


When  (Ks)max  is  very  large  relative  to  (Ks)mjn,  equation  (12)  is  approximately, 


K 


m 


=7  2  (K  ) 


min 


(14) 


For  the  maximum  fractional  power  change  over  the  range,  equation  (13)  becomes, 


K. 


y  max 


2  (K  ) 

v  s'  min 


g 


(15) 
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Figure  11  shows  the  table  with  several  ranges  of  thermal  conductivity;  the  pre¬ 
ferred  meter-bar  for  each  range,  and  the  maximum  fractional  power  change  associated 
for  a  given  geometrical  factor  of  1.0,  and  an  insulation  thermal  conductivity  of 


.058 


BTU 

HR-FT°F 


To  illustrate  the  fractional  power  change  that  would  arise  if  a  range  of  materials 
for  the  specimen  bar  were  measured  using  the  same  meter-bar  and  insulation  annulus, 
y  is  plotted  against  dimensionless  length  in  Figure  12.  The  geometrical  dimensions  of 
the  cut -bar  apparatus  in  Figure  12  are  the  dimensions  used  for  the  F  values  to  obtain 
a  plot  of  y  versus  the  dimensionless  length  u .  The  quantity  y  (u)  is  shown  for  specimens 
having, a  thermal  conductivity  ranging  from  one-tenth  to  ten  times  that  of  the  meters. 

The  following  is  an  example  problem,  illustrating  the  use  of  the  design  charts  and 
equations  of  a  cut -bar  apparatus  having  an  optimal  design  using  a  minimal  fractional 
power  change.  The  fractional  power  change  is  dependent  on  both  the  geometrical  factor 
and  the  thermal  conductivity  factor,  while  the  geometrical  and  thermal  conductivity 
factors  are  independent  of  one  another.  Thus,  these  factors  are  each  minimized 
separately. 

Assuming  that  the  particular  design  requirements  are  limited  to  an  overall  length  of 
w  equal  to  9.5  inches,  with  the  specimen  dimensions  of  length  L  equal  to  2  inches,  and 
specimen  radius  A  equal  to  1  inch,  the  specimen  thermal  conductivity  will  have  a  range 
from  a  minimum  of  4.5  BTU  to  a  maximun  of  242  BTU  (the  range 

of  most  metals).  HR-FTOF  HR-FT°F 


From  the  above  known  dimensions:  -jr-  =  2  and  =  .2111.  The  design 

**  B  W  TD 

chart  in  Figure  9  indicates  the  dimensionless  ratio  c  is  found  to  be  c  =  4.06 


or  Bc  =  4.06  inches.  Thus,  the  guard  diameter  should  be  at  least  4  inches  or  greater. 

The  computer  program  in  Appendix  B  is  used  to  calculate  the  geometrical  factor  for 
20  points  along  the  cut-bar.  The  computer  program  stops  when  the  geometrical  factor 
has  converged  within  the  test  range  of  .0001  or  reaches  the  maximum  number  of  series 
terms  of  M  equal  to  150.  The  frequency  F  of  the  test  for  convergency  occurs  every  7 
terms  in  the  series. 

The  inputs  into  the  computer  program  are:  W  =  9.5  inches;  A  =  1.0  inch; 

B  =  4.06  inches;  L  =  2.0  inches;  Test  =  .0001;  P  =  20  points;  M  =  150 
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terms,  and  F  =  7.  The  resulting  output  from  the  computer  program  is  illustrated  in  Figure 
13.  Calculated  for  each  point  is  the  distance  Z  of  that  point  in  inches  from  the  reference  end, 
the  number  of  KL  terms  before  the  last  term  in  the  series  is  calculated,  the  value  of  the 
SER  series  summation,  and  the  value  for  the  geometrical  factor  Fg.  #  The  resulting 
output  for  20  points  along  one  symmetrical  half  of  the  cut -bar  is  presented  in  Figure  13. 

From  this  output  of  the  computer  program,  a  graph  (Figure  14)  is  drawn  for  F  versus  Z  . 

In  Figure  14,  F  is  on  the  left  scale  of  the  ordinate  axis  and  its  maximum  value  at  the 
center  of  the  cut -bar  is  4.31. 

Since  the  apparatus  is  to  be  used  to  measure  the  thermal  conductivity  of  specimens 
ranging  from  Bismuth  to  that  of  Silver,  the  limiting  conductivity  values  are  (Kg)  .  =  4.5 


BTU 

HR-FT°F 


and  (Ks) 


s'  max 


=  242 


BTU 


HR-FT°F 


The  insulation  best  suited  for  the  appa¬ 


ratus  has  a  value  of  =  .  1 


BTU  .  Using  equation  (12),  Km  is  calculated  as  9 


HR-FT°F 


_ BTU  ,  which  is  in  the  range  of  stainless  steels.  Using  the  above  values  for  K., 

HR-FT°F  1 

(Ks)  min.  (Kg)  max  and  (Fg)  rnay  and  substituting  these  values  into  equation  (13),  the 

maximum  value  of  the  fractional  power  change  is  047  which  is  at  the  center  of  the 
specimen  bar. 

To  calculate  the  fractional  power  change  for  a  specific  test  the  apparent  con- 

RTT7 

ductivity  is  first  used,  for  instance  K  =  90 T~J~zTCr  .  Knowing  K.  and  K  the 
conductivity  factor  is  found  to  be  .01  from  equation  (8).  Since  the  geometrical 
factor  has  been  plotted  for  the  apparatus  in  Figure  14,  the  fractional  power  change 
can  be  found  by  changing  the  scale  on  the  axis  of  the  ordinate  by  a  factor  of  .01,  and 
this  scale  is  shown  on  the  right  scale  on  the  axis  of  the  ordinate. 

Since  the  temperature  gradients  are  calculated  between  two  thermocouples  placed 
a  finite  distance  apart,  tjie  fractional  power  changes  y  and  y  which  are  used  in 
equation  (4)  represent  the  average  fractional  power  changes  between  the  measuring 
stations,  this  average  y  can  be  calculated  as, 


y  aver 


y  dZ 


Z 


1 
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To  find  the  average  y  to  be  used  for  y  in  this  example,  (referring  to  Figure  14) 

the  cross-hatched  area  under  the  curve  from  Z^  =  1.1875  to  =  3.325  inches,  is 

measured  by  counting  graph  squares  with  the  result  that  y  =  .01025. 

y  s  is  calculated  by  using  the  same  graphical  procedure  with  the  two  measuring 

thermocouples  being  at  points  =  4.037  and  Z^  =  5.463  on  the  specimen  bar;  thus 

y  is  .  04 14 . 
x  s 

The  fractional  power  change  for  y  and  y  is  then  substituted  into  equation  (4), 

m  s 

BTU 

where  =  9  FT°F — %  The  resu^ing  equation,  to  find  the  thermal  conduc¬ 

tivity  of  the  specimen  for  the  example  problem,  is, 


S  and  S 
m  s 


Ks  =  8-72  r21 


s 

are  determined  from  the  measurements  of  the  thermocouples. 
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DESIGN  CHART 


DESIGN  CHART 

RATIO  B 
A 


FOR  TW  OPTIMUM  DIMENSIONLESS 
VARYING  THE  DIMENSIONLESS 


VARIABLES  L  AND  L 


FIGURE  10 
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EFFECT  OF  SPECIMEN  THERMAL  CONDUCTIVITY  ON  THE 

THERMAL  CONDUCTIVITY  FACTOR  FOR  A  RANGE  OF 
DIMENSIONLESS  RATIOS 


100  500  1000  500  10000 
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Maximum  Optimal  K 


Maximum  Error 
in  Kg  Range 
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BTU 


HR-FT°F 


0.3872 


9.0753 


10.5183 


11.4430 


105.1839 


FIGURE  12 


PLOT  OF  FRACTIONAL  POWER  CHANGE  FOR  A  RANGE  OF  SPECIMEN 
THERMAL  CONDUCTIVITIES  VERSUS  DIMENSIONLESS  LENGTH 


Dimensionless 
apparatus  are:  i 


parameters  for  this 
*2;  |  =  .2105;  B  = 


3 
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FIGURE  13 
COMPUTER  OUTPUT 


RESULTS  OF  POWER  LOSS  OVER  LENGTH  OF  METER- BAR 


F  VERSUS  Z  FOR  20  POINTS  FROM  Z=0  TO  END 
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PLOT  OF  THE  GEOMETRICAL  FACTOR 


THE  LEFT  SCALE  AND 


FRACTIONAL  POWER  CHANGE  ON  IKE  RIGHT  SCALE  OF  THE 
ORDINATE  AXIS  VERSUS  THE  Z  ON  THE  AXIS  OF  THE  ABSCISSA 

Dimensions  of  this  apparatus  are; 

W  —  P 1 5  in  e «  A  I  in  on  B  ®  4o  06  m  no  L  2  b  in  o « 
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APPENDIX  A 


MATHEMATICAL  ANALYSIS  OF  GEOMETRICAL  FACTOR 


It  is  important  to  note  that  this  analysis  is 
based  on  the  change  of  heat  flow  in  the  annulus  and 
through  the  boundaries  of  the  annulus,  so  the  equa¬ 
tions  and  solution  are  based  only  on  the  annulus,  and 
not  on  the  cut-bar  or  the  cylindrical  guard; 
however,  the  boundary  conditions  are  based  on  the 
temperature  distribution  of  the  cut-bar,  the  cylindri¬ 
cal  guard,  and  the  heater  cooler  system. 

Mathematical  investigation  of  the  powder  insula¬ 
tion  annulus  follows: 

Figure  A-l  shows  a  schematic  of  the  cut-bar 
apparatus.  The  two  meter  bars  are  of  equal  length 
M  and  of  the  same  thermal  conductivity  K^. 

W  -  the  overall  length  of  the  cut-bar 

M  -  the  length  of  the  meter-bar 

L  -  the  length  of  the  specimen  bar 

A  -  the  radius  of  the  cut-bar,  or  the  inner 
radius  of  the  powder  insulation  annulus 

B  -  the  outer  radius  of  the  powder  insulation 
annulus  or  the  inner  radius  of  the  cylin¬ 
drical  guard. 

K  -  the  thermal  conductivity  of  the  meter-bar 
m  y 


26 


la! 


A  SCHEMATIC  OF  THE  CUT- BAR  APPAR* 


Kg-  the  thermal  conductivity  of  the  specimen  bar 

K.-  the  thermal  conductivity  of  the  powder  insu¬ 
lation  in  the  annulus 

V  -  the  temperature  variable 

The  boundary  conditions  of  the  powder  insulation 
annulus  are: 


r  =  A 

0  i  z  ^  w 

V  =  G  (Z) 

r  =  B 

0  *  Z  S  W 

V  =  H  (Z) 

A  =  r  =  B 

Z  =  0 

V^=  V  =  Constant 

G(0)  =  H(0)  =  VL 

PQ 

Nl 

U 

Ml 

< 

z  =  w 

V  =  V2=Constant 

G(W)  =  H(W)  =  V2  (A-l) 

The  above  boundary  conditions  of  the  annulus  may 
be  reduced  to  dimensionless  variables. 

The  dimensionless  variables  are: 


u  =  Z 

m  =  M 

»  =  L 

W 

W 

W 

P  =  r 

a  =  A 

b=  B 

w, 

W 

W 

e  =  v  -v2 

_  K 

_  K. 

(A-2) 

vl-v2 

Km 
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The  dimensionless  boundary  conditions  are: 
/=  a  0  -  u  -  1  e  =  g(u)  =  G(Z)~V2 


V1  -  V2 


J>=  b 


0  -  u  -  1  9  =  h(u)  =  H(z)~^2 


V1  -  V2 


a  -  P  ~  b  u  =  0 

a  =  f  =  b  u  =  1 


0  =  1 
0=0 


(A- 3  ) 


The  governing  differential  equation  of  heat  flow 
in  the  annulus  is  Laplace* s  equation  in  cylindrical  co¬ 
ordinates  : 

v2v  =  o 

Using  dimensionless  variables  Laplace’s  equation 
is. 


S>20  l  90  S>a© 

ap5  p  9P  +  3  u?-  U 

(A-4) 

Using  separation  of  variables, 

e  =  pe 

(A-5) 

Substituting  equations  (A-5)  into  (A-4) : 

p.« 


Therefore, 


And , 


%  +  k  £'  +  S"  -  o 

p  y  p  e 


E"  =  - 

E 


p"  +  1  pi  .  =  0 

P 


(A-6) 

(A-7) 
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Solving  equation  (A-6), 

E  =  C^Sin(«<u)  +  C2Cos(<*u)  (A-8) 

Solving  equation  (A-7)  the  solution  is  in  the 
form  of  Modified  Bessel  Functions: 

(A-9) 


P  =  C3I0(eip)  +  C4Kq(o<P) 
Substituting  equations  (A-8)  and  (A-9)  into 


(A-5) : 


9  = 


C3I0(^P)  +  C4K0(«*P)  |Lsin(*u)  +  C2Cos(^u) 

-  —II—  Ta^io) 


Applying  the  last  boundary  condition  to  (A-10): 
a=P=b  u  =  1  0=0 


Then, 


=  0  Sin<*=  0  and  <A  =  nTT 


Where,  n  =  0,  1,  2,...  thus,  equation  (A-10) 
becomes : 


00 

6  -  C0  +5=1 


— I 


CnI0(n1tf’)+DnK0(n'tTS>)  Sin(nlTu)  (A-ll) 


Using  the  other  boundary  conditions  from  equations 


(A-3): 


a  -  ?  -  b 

c 

ii 

o 

(D 

II 

h-1 

a  =  ?  =  b 

u  =  1 

O 

ii 

<D 

Then, 

CQ  =  1-u 
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And , 


oo 

0=l-u+  nS  L  [cnI0(n'n>»)+DnK0(n'rrp) 


Sin(nTTu)  (A-12) 


At  the  inner  diameter  of  the  annulus,  using 

equation  (A-12)  the  equation  for  the  inner  perimeter 

surface  temperature  is: 

00 

g(u)=  (l-u)+E  1  |cnI0(nTTa)+DnK0(n'tra  )  J  Sin(nUu) 


(A-13) 


Correspondingly  at  the  outer  diameter  of  the  annu¬ 
lus  using  equation  (A-12)  the  equation  of  the  outer 
perimeter  surface  is: 

00  -j 

h(u)=  (l-u)-^S  ^  |  CnlQ(nTrb)+DnKo(nTr  b)  J  Sin(ntr  a) 

( A-14) 


Letting, 


an  =  CnIQ  (n  TV  a)  +  DnKQ  (n*ra) 


(A-15) 


bn  =  Vo  (nTTb)  +  DnK0  <n*b> 


(A-16) 


Solving  C  and  D  from  equations  (A-15)  and  (A-16): 


n 


n 


C  = 
n 


an 


n 


*0  (n  fra) 


0  (n  rrb) 


K. 


T0  Cnrr  a)  Ko  (ntTa) 

X0  (ntTb)  *0  (nTTb) 


(A-17) 
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I0(nTTa)  an 

In(nirb)  b 

u  n 

IpCnTTa)  K^nira) 

I0(nttb)  K0(nT(b) 


(A-18) 


Substituting  equations  (A-15)  into  (A-13)  and 

(A-16)  into  (A-14)  and  arranging  so  that: 

00 

u-l+g(u)=  T*  a^SinCnTTu)  (A-19) 

n=L  n 
00 

u-l+h( u)=  ^  bnSin(nTTu)  (A-20) 


Solving  for  the  Euler  coefficients  a^  and  bn,  it 
is  necessary  to  use  Euler-Four ier  formulas,  ^ 


The  ha If -range  sine  expansion  is  based  upon  ex¬ 


tending  (u-l)+g(u)  and  (u-l)+h(u)  over  the  interval  u 


of  (-1,0)  by  reflection  in  the  origin. 


1 :  Drawing 


If  BB^  is  chosen  as  the  vertical  axis,  the  graph 
defines  an  odd  function  by  the  theorem  only  sine  terms 
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If, 

p  =  1  and  u=  t  (from  footnote  1) 

a  and  b  are : 
n  n 


an  =  2  J  ^(u-l)+g(u)j  Sin(nTtu)du 
bn  =  ^ Jq  ^(u-l)+h(u)J  Sin(nTru)du 


(  A-21) 
(A-22) 


Using  the  boundary  conditions 


for  a  linear  guard: 


h(z)=v1-(vl-v2)z 

Substituting  into  equation  (A-23)  the  dimen¬ 
sionless  variables  of  the  equation  (A-2)  and  the 
boundary  condition  h(u)  of  equations  (A-3): 

h(u>  = 


(A-23) 


V1'V2 


Thus , 


h(  u)  =  1— u 

By  substituting  equation  (A-24)  into  (A-22): 


bn  =  0 


(A-24) 


(A-25) 


will  appear  in  its  expansion. 

Theorem:  If  f(t)  is  an  odd  periodic  function,  then  the 

coefficients  in  the  Fourier  series  f(t)  are  given  by 
the  formulas: 

A  =  0  B  =  2  f  P  f  ( t)Sin(nTT t)dt 

n  n  p  Jq  p 

The  oddness  depends  upon  the  relation  to  the  vertical 
axis  of  the  coordinate  system  which  can  be  arbitrarily 
chosen. 
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Differentiating  equation  (A-12): 

d6  =  V  nTTSin(rfnu)  [ CnIo(rOTP)+DnKo(n'nf!)"| 
df  n=l  L 

=  XI  nTTS in(nllu)  |^CnI1(n1TP)-DnK1(nlT?)J  (A-26) 
n=l 


Let  equation  (A-26)  define  the  surface  of  the  annu¬ 
lus  at  the  inner  diameter  A. 

Substituting  equations  (A-17) , (A-18) , (A-25) 

into  (A-26)  the  resulting  equation  for  the  inner 

diameter  is :  _ 

K0(nTTb)I^(njfTa)+I0(nn'b)K^(n'nra) 

KQ  (ntrb)  IQ  (ntta )  -i0(nirb jKQ(rmaj 

(A-27) 

The  radial  heat  flow  at  the  inner  surface  of  the 
annulus,  r  =  A,  through  a  cylindrical  surface  ele¬ 
ment  of  length  dZ  is: 


oo 

4f=XXnffSin(mru) 

dJr  n  =  1 


dp  =  2TT  A  K^SV  j 


dZ 


r=A 


(A- 28) 


Using  the  dimensionless  variables  from  the 
equations  (A-2) ,  equation  (A-28)  becomes: 

dp  =  2TTac-:  /33\  du  (A-29) 

The  total  dimensionless  radial  heat  flow  across 
the  inner  cylindrical  surface  of  the  annulus  is 
found  by  integrating  equation  (A-29)  for  limits 
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between  u  =  to  u  =  u>  which  are  defined  in  the 
analysis  subsequently. 


(A-30) 


Substituting  equation  (A-27)  into  (A-30)  and  integrating 
oo  r  1 

p(u)^  2TTa  (r^an|cos(nTTuI)-Cos(nTru2)J 

K0(nrT  b)I^(nlTa)+lQ(nTT  b)K^(nlTa) 
l(0CnTr  b;i(-)(ln'n'a;-I0(niT  b;K!0Cn7ra^ 

(A-31) 

It  is  necessary  to  define  a  temperature  distribu¬ 
tion  along  the  inner  surface  of  the  annulus,  which 

then  defines  a  . 

n 

Referring  to  Figure  A-L:  Vi  is  the  temperature  at 

the  matched  end  where  Z  =  0;  S  is  the  constant  longi- 

s 

tudinal  temperature  gradient  in  the  specimen  bar;  and 
is  the  constant  longitudinal  temperature  gradient 
in  the  meter-bar. 

The  boundary  conditions  of  the  cut-bar  define  the 
temperature  distribution  along  the  inner  surface  of 
the  annulus.  This  temperature  is  constant  around  the 
perimeter  which  is  equal  distant  from  the  ends. 
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0  -  Z  -  M 
M  -  Z  •  M+L 

m  +  L  -  Z  -  W 


V  =  T-SZ 

m 

V  =  T-S  M-S  (Z-M) 

HI  s 

V  =  T-S  M-S  L-S  (Z-M-L)=S  (W-Z) 

m  s  m 

(A-32) 


Using  the  dimensionless  variables  of  equations 
(A-2)  the  above  boundary  conditions  (A-32)  are  re¬ 
defined  for  a  dimensionless  temperature  distribution 
g(u)  : 

O-u-m  0  =  1-  ¥  u 

m 


m  -  u  -  m  +  8 

e 

=  1-  H*  m-  V  (u-m) 

m  sN 

m  +  8  -  u  -  1 

e 

=  S’  <l-u) 
m 

(A-33) 

Where , 

s  w 

S  W 

Y  = 

m 

,  »=  s 

(A-34) 

1  m 

T 

¥  x 

Using  the  above  boundary  conditions  (A-33)  and 
equation  (A-21)  the  value  for  can  be  defined  as: 


anl= 

anl= 


an  ~  anl+an2+an3 


if*  u(  1  —  4»m)Sin(n1Tu)du 

2(  1-  \,)|sin(n  rt  m)  -n'rrBCos(n«m)l 
n2T T2 


(A-35) 


(A-36) 
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I+S 


a  _  =  2 |  (-  ¥  m+  T  m)Sin(nTTu)du 

n2  J  ms 

•'m 

+2  J  (l-  )uSin(nTTu)du 
s 

a  o  ^  %n)  [" Cos(irtfm)-Cosntr  (m-  Jl  )1 

n  nTT  L  J 

2(1-  4*  )  r 

+  - * — is-  Sin(n  rr  m+nTT^ )-n  TT(nH- J?  )Cos  (ntTm+rrrr^) 

ri  1TZ  L 

-Sin(nrrm)+n  TTrnCosCntTmjj  (A- 37 ) 

an3‘2  fl  [V1*"'1- V]  Si''<"mJ,du 
*/B  +  S 

an3~2(Ym~ J  )  j  Cos  [  nrr  C rrH-JJ)]  -Cos(nTr)  | 


2(1-  )  r 

+•  —n  ■  TO  -n'CTCosCnTT  )-Sin(n1T  nrl'nfrJ!) 
Tr  *  L 

+ntr(m+ft)  Cos(nTTm+n*fT.4  )J 


(A-38) 


Substitute  equations  (A-36),  (A-37),  (A-38) 
into  (A-35)  and  by  arranging  and  reducing  the  Euler 
coefficients  become: 


an  "  — 0--a  (^m  -  H's)  [  Sin(--^  ^  Cos  (nTT )  1  (A-39) 

nlTC  L  Z  2  J 

From  equations  (A-2)  of  dimensionless  variables, 

(A-40) 


cr  =^s 

S  Km 
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Assuming  that  there  is  no  heat  loss  or  gain  in 
the  cut-bar,  it  is  possible  to  help  reduce  the 
mathematical  analysis  so  that  the  resulting  solution 
can  be  solved  by  computer  with  only  a  minor  loss  in 
the  accuracy  of  the  resulting  data. 


Thus , 

K  S  =  K  S 
mm  s  s 


Or, 


( A-41) 


From  equations  (A-34),  and  (A-41)  it  can  be  shown 
that , 


Vs  = 

KV  K 

Tm  s 


(A-42) 


Let  Q  represent  the  total  Longitudinal  heat  flow 
through  the  cut-bar  with  the  stipulation  that  there 
be  no  heat  gained  or  lost  through  the  surface  between 
the  cut-bar  and  the  annulus. 

Thus , 

Q  =  TTA2K  S  =  tTA^  S  (A-43) 

x  mm  s  s 


Defining  heat  flow  Q  in  dimensionless  parameters: 


(A-44) 
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Using  equations  (A-34),  (A-43)  and  (A-44),  we 


have , 


IT  V 


ra 


(A-45) 


Defining  the  net  fraction  of  power  lost  or  gained 
in  the  annulus  between  u  =  0  and  u  =  Uj  as  fl  ,  Heat 
loss  in  the  cut-bar  is  the  heat  gained  in  the  annulus; 
thus,  the  fractional  power  change  V  in  the  cut-bar  is, 

ts  =  -9  =  -  |  (A-46) 

Substituting  in  values  for  p,  q  into  equation 
(A-46)  by  using  equations  (A-31) ,  (A-39),  (A-42), 
and  (A-45): 

*  =  -1 


/l_  -  l_\_j 

X 

l-Cos(nTT  u  ^ 
T 

Ik  k  L; 

'  m  s  'TV' 

a  £i 

n2 

K0(nTT  b)l^(n  fr  a)+I0(nTrb)K^(nTr  a) 

Sin/n'tTJ?  ^  Cos  /  nfr^ 

K^(ntr  b)I^(nTT  a) -I^(nnb)K^(nTf  a) 

V  2  )  {  2  ) 

(A-47) 

Substituting  the  dimensionless  variables  of 
equation  (A-2)  into  equation  (A-47),  equation 
(A-47 )becomes : 


00 


Jf  =  -K. 


(l_  -  l\  _8W  V 

AKS  KsJ  7T2A 

n=l 


1-Cos 


n 


KC 

Ki( 

p  nrr  A\ 

.  w  / 

KC 

fnirA\_ 

T 

Ko< 

TTA^ 

s,  w  / 

0\  W  J 

^  w  / 

\  w  / 

(A-48) 
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X  is  the  per  cent  heat  loss  relative  to  the 
total  heat  flow  in  the  cut-bar  assuming  that  the 
total  heat  flow  is  the  flow  when  there  would  be  no 
loss  or  gain  in  the  flow  from  or  to  the  powder 
insulation  annulus* 


%  is  called  the  fractional  power  change. 

Letting, 

#  =  F.F 
k  g 

Where, 


(A-49) 


(A-50) 


F^is  the  thermal  conductivity  factor  of  the 
fractional  power  change, 

And, 


oc 


F  =  8£L 

g  TT^Ai 


1-Cos  ( 

nTTz\ 
w  / 

I  2 

i  n 

v  fSJEL§ 
KcA  W 


v  /nrr  B 

ot“w“ 


)h¥v± 


)♦ 


Sin 


(#)  c°*(rt} 

(A-51) 


It  is  necessary  to  simplify  equation  (A-51) 
for  computer  analysis  since, 

Cos(n^  }  -  0  when  n  is  odd 
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Therefore, 


n 

(-1)^  =  Cos(n-^  )  when  n  is  even  (A- 52) 

Letting,  m  =  ^  (A-53) 

Substituting  equations  (A-53),  (A-52)  into 
equation  (A-51),  the  index  of  equation  (A-5J.)  can 
be  changed  with  no  Loss  to  the  equation. 

Therefore, 


00 


F  is  the  geometrical  factor  of  the  fractional 
S 

change. 
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APPENDIX  B 

COMPUTER  PROGRAM  OF  GEOMETRICAL  FACTOR 


A  computer  program  was  developed  to  calculate  the 
large  number  of  data  points  needed  to  plot  the  design 
charts.  Also,  this  program  was  utilized  in  the  error 
analysis  of  any  particular  apparatus,  A  particular 
range  of  design  dimensions  can  be  identified  by  the 
use  of  these  charts.  The  exact  error  for  a  given  set 
of  dimensions  can  be  found  by  using  the  program. 

In  Appendix  A,  the  resulting  solution  is  described 
in  a  mathematical  equation  (A-54)  for  the  geometrical 
factor  of  the  cut-bar  apparatus.  This  equation  is 


reproduced  below: 

00 


=  2w  y\_ 

'  ti2a  ' 


1) 


m 


m  =  1 , 2  ,  . , 


1-Cos 


m 


*I0  (laTI°)K1(aSffiA) 
-I0  <^P>K0<^> 


Sin  (SOIL) 
W 


(B-l) 


Two  computer  programming  approaches  are  obvious. 
The  first  approach  uses  a  large  memory  bank  where  all 
of  the  Modified  Bessel  Function  of  the  first  and 
second  kind,  and  of  the  first  and  second  order,  are 
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stored.  Using  this  method,  it  is  necessary  to  use 
approximation  techniques  to  find  the  Modified  Bessel 
Functions  of  fractional  arguments.  This  method  can 
be  used  only  in  a  data  reduction  center  having  the 
necessary  facilities.  The  second  approach  is  to 
generate  the  Modified  Bessel  Function  in  a  program 
when  needed. 

The  second  approach  was  used  in  this  thesis  be¬ 
cause  the  IBM  1620  Computer  which  was  available  was 
limited  in  size  and  speed,  and  the  computer  program 
using  the  Fortran  language  was  efficient.  The  pro¬ 
gram  can  function  in  the  Fortran  II,  Fortran  IV,  or 
Forgo  language,  and  with  a  small  change  it  can 
function  in  PDQ  language. 

Inputs  into  the  following  programs  are:  W,  A, 

B,  CL,  TEST,  JA,  MM,  LP. 

W  -  length  of  the  annulus. 

A  -  the  diameter  of  the  specimen  bar  and  meter- 
bar,  or  the  inner  diameter  of  the  annulus. 

B  -  the  outer  diameter  of  the  annulus  or  the 
inner  diameter  of  the  cylindrical  guard. 

CL-  the  longitudinal  length  of  the  specimen 
bar.  It  is  represented  by  L  in  the 
equation  (B-l). 

TEST-  the  accuracy  of  the  resulting  output 

TEST  stops  the  program  or  truncates  the 
series  when  it  has  converged  so  the  Fg  is 

varying  within  the  range  of  the  TEST  input. 
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JA  -  the  number  of  points  for  which  Fg  output 
is  formed  along  one-half  of  the  bar.  It 
is  the  number  of  Z  points  into  which  one- 
half  of  the  bar  is  dividided. 

MM  -  the  limiting  number  of  terms  in  the  series. 
It  represents  the  maximum  number  of  terms 
at  which  the  series  or  program  will  be 
truncated.  At  some  points,  the  series  will 
converge  only  at  a  very  large  number  of 
terms . 

LP  -  the  frequency  of  the  test  for  convergence 
of  the  series.  It  must  be  an  odd  number 
and  it  is  based  on  the  2W.  with  familiarity 

L 

in  using  the  program,  an  engineer^  intui¬ 
tiveness  can  find  a  reasonable  LP. 


Figure  B-l  illustrates  the  flow-chart  for  the 
computer  program.  Three  computer  programs  were 
developed  based  on  this  flow-chart  and  the  available 
facilities.  The  computer  program  illustrated  in 
Figure  B-2  was  developed  for  use  on  an  IBM-1620  com¬ 
puter.  Using  this  program  as  a  source  deck,  an 
objective  deck  can  be  compiled  in  Fortran  II  and  FDQ. 
This  program  can  be  used  directly  in  Forgo  language. 
Figure  E>-3  illustrates  a  computer  program  which  can 
be  used  on  a  t ime-ehar ing-computer  system,  using 
Fortran  IV  language,  and  this  program  can  be  modi¬ 
fied  with  the  addition  of  the  program  statement  lines 
illustrated  in  Figure  B-4.  This  modification  expands 
the  program  so  that  a  large  amount  of  data  can  be 
hand  led . 
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FIGURE  B-1 


COMPUTER  PROGRAM  FLOW-CHART 


FIGURE  B-2 


COMPUTER  PROGRAM  FOR  IBM- 1620  COMPUTER 

DIMENSION  FAB( 20) ,FB(20) ,Z(20) ,FG( 20) , SER( 20) , 
BST( 20) ,KL( 20) 

PI=3. 1415927 

GAM=. 57721566 

PUNCH  910 

10  READ  905,W,A,B,CL,TEST,JA,MM,LP 
PUNCH  920, JA 
PUNCH  930 

PUNCH  935 ,W, A, B,CL,TEST, JA,MM,LP 
SIGN  =  1, 

EP  =  1. 

DOG  =  1. 

ENT  =  EXPF((-4.*(B-A)*PI)/W) 

H  =  ( 2 . 0*W) / ( PI*PI*A) 

CAT  =  JA 
DO  20  N=  1 , JA 
ZA  =  N 

Z  (N)  =  ZA*W/(2.*CAT) 

SER(N)=0 , 0 
BST(N)=9999. 

KL(N)=0 

20  FAB(N)=1 . 0 
M=0 
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30  M=M+1 

DOG=DOG*ENT 

V=M 

XV=2.*V*PI/W 
RRI=V*PI/CAT 
SIGN=  -SIGN 

Y=SINF(V*PI*CL/W)*SIGN/(V*V) 

X=A*XV 

CKS=X 

KK=0 

40  NORD=0 

500  BK=0. 0 
FN=NORD 

IF(X-FN-6 . ) 501 , 700 , 700 

501  XA=X/2  . 

XB=XA*XA 

N=0 

503  AN=N 
T=1 . 

S=  -1. 

510  IF (AN) 9999 ,520,512 
512  T=T*XA/AN 
AN=AN-1. 

GO  TO  510' 

47 


520  BIN=T 

DO  530  K=l,9999 
DEN=K*(K+N) 

T=T*XB/DEN 

IF ( ( BIN+T) -BIN)525 ,550,525 

525  BIN=BIN+T 

530  CONTINUE 

550  IF(N-1)575,630,555 

555  IF(X-FN-3.) 1111, 1111 , 700 

565  N=1 

BIO=BIN 

BKO=BK 

GO  TO  503 

575  BK=  - ( GAM+LOGF (XA) )*BIN 
T=XB 
8=1. 

XI=1. 

DO  610  K=2 , 9999 
AK=K 

IF((BK+T*XI)-BK)600 ,620,600 
600  BK=BK+T*XI 

T=T*XB/(AK*AK) 

XI=XI+1./AK 
610  CONTINUE 
620  IF(NORD)9999 , 555 , 565 
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630  BK= ( 1 . /X- BIN* BKO ) / BIO 
IF(NORD- 1)999 6, 555 ,9999 
700  C=4*N0RD*N0RD 
D=8.*X 

C0N2=1./SQRTF(2.*PI*X) 
C0N3=SQRTF ( PI / ( 2 . *X) ) 
AN=NORD 

PHI=X- ( 2 . *AN+1 . ) /4 . *PI 
K=X+ 1 . +SQRTF ( X*X+AN*AN ) 
T=(C-1.)/D 
S=l. 

U=l. 

PN=1. 

QN=T 

BK=1 . +T 

BI=1.-T 

DO  735  1=2, K 

AI=I 

T=(C-(2.*AI-1.)**2)/D*T/AI 

BK=BK+T 

BI=BI+T*S 

IF((BI+T)-BI >730,738,730 
730  S=  -S 
735  CONTINUE 
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738  IF(CKS-7. >741,755, 755 
741  BK»£XPF(-X)*C0N3*BK 

IF(X-FN-6.) 1111, 750, 750 
750  BIN=EXPF(X)*CON2*BI 
GO  TO  1111 
755  BK=C0N3*BK 
760  B1N=C0N2*B1 

EP=DOG 

1111  N0RD=FN 

IF (NORD) 9999, 810, 830 
810  IF(KK)9999, 820,840 

820  GA=BIN 
GB=BK 
NORD=l 
GO  TO  500 
830  GC=BIN 
GD=BK 
X=XV*B 
KK=1 

GO  TO  40 
840  GE=BIN 

gf=bk*ep 

R=( GF*GC+GE*GD ) / (GF*GA-GE*GB) 
FA=Y*R 
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LPF=0 


DO  860  N=l, JA 
IF (KL(N) -3 >845,850,858 
845  LPF=1 
850  ZA=N 

FB(N) =1 . O-COSF (RRI*ZA> 

FAB(N)=FA*FB(N) 

ser(n)=ser(n)+fab(n) 

IF (M/LP*LP-M) 860 , 852 , 860 
852  IF(TEST-ABS( BST(N)-SER(N) ) >854,856,856 
854  KL(N>=0 

BST(N)=SER(N) 

GO  TO  860 
856  BST(N)=SER(N) 

858  KL(N)=KL(N>+1 

860  CONTINUE 

IF(LPF-l) 890,888,9999 
888  IF(M-HM)30, 890 ,890 

890  PUNCH  960, M 
DO  900  J=1 , JA 
FG( J>=H*SER( J) 

900  PUNCH  940, J,Z  ( J) ,KL( J) ,SER( J) ,FG( J) 

GO  TO  10 

9999  PUNCH  9998 
GO  TO  10 

905  FORMAT  (5F10.5,3I5) 
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910  FORMAT  ( ///12X,46HRESULTS  OF  POWER  LOSS  OVER 
LENGTH  OF  METER- BAR//) 

920  FORMAT  ( ///12X, 15HFG  VERSUS  Z  FOR,l3,23H 
POINTS  FROM  Z=0  TO  END//) 

930  FORMAT  (6X. 1HW.9X, 1HA, 9X, 1HB,9X, 1HL, 8X.4HTEST, 
5X , 1HP , 3X , LHM, 4X , 1HF ) 

935  FORMAT  ( 5*10 . 5 , 315/ // ) 

940  FORMAT  (110 ,F14. 5 , 18 , 3X, 2F10 . 7) 

960  FORMAT  (//7X,5HP0INT,8X,1HZ,9X,2HKL,7X,3HSER, 
8X, 2HFG , 5X, 2HM= , 13) 

9998  FORMAT  ( //25X, 5HERR0R//) 

END 


52 


FIGURE  B-3 


COMPUTER  PROGRAM  FOR  TIME- SHARING  COMPUTER 

10  DIMENSION  FAB(30),FB(30),Z(30),FG(30) ,SER(30) , 
BST(30) ,KL(30) 

20  PI=3. 1415927 

30  GAM=. 57721566 

40  PRINT  910 

50  10  INPUT,W,A,B,CL,TEST,JA,MM,LP 

60  PRINT  920, JA 

70  PRINT  930 

80  PRINT  935 ,W,A, B,CL,TEST, JA,MM,LP 

85  SIGN=1 . 

90  EP=1. 

100  DOG=l. 

110  ENT=EXPF((-4.*(B-A)*PI)/W) 

120  H=(2.0*W)/(PI*PI*A) 

130  CAT=JA 

140  DO  20  N«1,JA 

150  ZA=N 

160  Z(N)=ZA*W/(2»*CAT) 

170  SER(N)=0, 0 

180  BST(N)=9999. 

190  KL(N)=0 

200  20  FAB(N)=1.0 
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210  M=0 

220  30  M=M+1 

230  DOG=DOG*ENT 
240  V=M 

245  XV=2.*V*PI/W 
248  SIGN=  -SIGN 

250  Y=SINF(V*PI*CL/W)*SIGN/(V*V) 

260  RRI=V*PI /CAT 

280  X=A*XV 

290  CKS-X 

300  KK=0 

330  40  N0RD=0 

340  500  BK=0 . 0 

350  FN=NORD 

360  IF(X-FN-6. >501,700,700 
370  501  XA=X/2. 

380  XB=XA*XA 
390  N=0 

400  503  AN=N 

405  T=l. 

410  S=  -1. 

415  510  IF (AN) 9999, 520, 512 

420  512  T=T*XA/ AN 

425  AN=AN-1 . 

430  GO  TO  510 
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435  520  BIN=T 

440  DO  530  K=L, 9999 

445  DEN=K*(K+N) 

450  T=T*XB/DEN 

455  IF( (BIN+T)-BIN)525i 530, 525 

460  525  BIN=BIN+T 

465  530  CONTINUE 

470  550  IF (N- 1)57 5, 630, 555 

475  555  IF(X-FN-3. >1111, 1111, 700 

480  565  N=1 

485  BIO=BIN 

490  BKO=BK 

495  GO  TO  503 

500  575  BK=  - ( GAM+LOGF ( XA) ) *BIN 

505  T=XB 

510  S-l. 

515  XI=1. 

520  DO  610  K=2, 9999 

525  AK=K 

530  IF((BK+T*XI)-BK)600,620,600 

535  600  BK=BK+T*XI 

540  T=T*XB/(AK*AK) 

545  XI=XI+1./AK 

550  610  CONTINUE 

555  620  IF (NORD)9999 ,555, 565 
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560  630  BK=(L./X-BIN*BKO)/BIO 

565  IF (N0RD- 1)9999, 555, 9999 
570  700  C=4*N0RD*N0RD 

575  D=8.*X 

580  CON2=l./SQRTF(2.*PI*X) 

585  C0N3=SQRTF( PI/ ( 2 . *X) ) 

590  AN=NORD 

600  PHI=X-(2.*AN+1.)/4.*PI 

602  K=X+ 1 . +SQRTF ( X*X+AN*AN ) 

604  T=(C-1 . )/D 

606  S=l. 

608  U=l. 

610  PN=1. 

612  QN=T 
614  BK=1.+T 

616  BI=1.-T 

618  DO  735  1=2, K 
620  AI=I 

622  T=(C-(2.*AI-1.)**2)/D*T/AI 

624  BK=BK+T 
626  BI=BI+T*S 

628  IF( (BI+T)-BI)730, 738, 730 
630  730  S=  -S 

632  735  CONTINUE 
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634  738  IF(CKS-7. >741,755,755 

636  74L  BK=EXPF(-X)*C0N3*BK 

638  IF(X-FN-6. >1111,750,750 

640  750  BIN=EXPF ( X ) *CON2*BI 

642  GO  TO  1111 

644  755  BK=C0N3*BK 

646  760  BIN  =€ON2*BI 

648  EP=D0G 

650  1111  N0RD=FN 

652  IF(N0RD ) 9999 , 810 , 830 

654  810  IF (KK) 9999 ,820,840 

656  820  GA=BIN 

658  GB=BK 

660  N0RD=1 

662  GO  TO  500 

664  830  GC=BIN 

666  GD=BK 

668  X=XV*B 

670  KK=1 

672  GO  TO  40 

674  840  GE=BIN 

676  GF=BK*EP 

678  R= ( GF*GC+GE*GD )/ ( GF*GA- GE*GB ) 

680  FA=Y*R 
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682  LPF=0 
1  684  DO  860  N=1,JA 
686  IF(KL(N)-3)845,850,858 

688  845  LPF=1 

690  850  ZA=N 

692  FB(N)=1.0-COSF(RRI*ZA) 

694  FAB(N)=FA*FB(N) 

696  SER(N) =SER(N)+FAB(N) 

698  IF(M/LP*LP-M)860,852, 860 

700  852  IF (TEST- ABS ( BST(N )- SER (N) )) 854, 856,856 

702  854  KL(N) =0 

704  BST(N)=SER(N) 

706  GO  TO  860 

708  856  BST(N)=SER(N) 

710  858  KL(N)=KL(N)+1 

712  860  CONTINUE 

714  IF(LPF-1)890, 888, 9999 

716  888  IF(M-MM) 30 , 890 , 890 

718  890  PRINT  960, H 

720  DO  900  J=1 , JA 

722  FG ( J ) =H*SER ( J ) 

724  900  PRINT  940, J,Z  ( J) ,KL( J) ,SER( J) ,FG( J) 

726  GO  TO  10 

728  9999  PRINT  9998 
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730 


GO  TO  10 


732  910  FORMAT  (//12X  46HRESULTS  OF  POWER  LOSS 

OVER  LENGTH  OF  METER-BAR//) 

734  920  FORMAT  (///12X.15HFG  VERSUS  Z  FOR,  13, 
23H  POINTS  FROM  Z=0  TO  END//) 

735  930  FORMAT  (6X,1HW,9X,1HA,9X,1HB,9X,1HL,8X, 

4HTEST , 5X , 1HP , 3X , 1HM, 4X , 1HF ) 

738  935  FORMAT  (5F10.5.3I5///) 

740  940  FORMAT  (U0,F14.5,18,3X,2F10.7) 

742  960  FORMAT  ( //7X, 5HPOINT, 8X, 1HZ, 9X,2HKL, 7X, 

3HSER , 8X , 2HFG , 5X , 2HM= , 13 ) 

744  9998  FORMAT  ( //25X, 5HERROR// ) 

746  END 
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FIGURE  B-4 

COMPUTER  PROGRAM  FOR  TIME-SHARING  COMPUTER 

50  10  READ ,W, A, B,CL,TEST , JA,MM,LP 
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